Swiss Glaciers

#rm(list = ls())

1 Introduction

Glaciers around continental Europe have been receding for some time now, and as the most mountainous country in Europe, Switzerland is acutely affected by this phenomenon. The warming climate has meant a change in the alpine nature in the country, as well as what this means going forward. Switzerland is famed for its nature, and its Glaciers are the crowning glory among its magnificent alps. This project aims to analyze glacial mass differences from 1850 to 2010, one hundred and sixty years of change, to highlight how much of the glaciers has melted, as well as the regions most affected. But one thing is clear, if you want to see glaciers, start making plans sooner rather than later.

2 Data Sources

In order to map out the areas that have been affected by Glacial Regression, data has been taken from GLAMOS, a publicly funded entity operated directly by the University of Zurich and Fribourg as well as ETH Zurich. Official government websites, as well as GADM for country specific shape files. The data gathered spans one hundred and sixty years from 1850 to 2010 showing glacial changes throughout multiple centuries. ,

3 Glacier Dense Regions

For better visualization, these are some drone videos and photos from parts of the Swiss Alps

4 Calculating Glacial Levels in 1850 and 2010

Understanding the Glacial Regression means isolating the two years furthest apart in the available data for the difference to be most clearly seen, this case it would mean comparing 1850 to 2010 and to see how much has changed.

4.1 Reading the Data

#Setting Directory
library("dplyr")

Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

library(sf)

Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUE

setwd("/Users/gabriel_caffrey/Documents/FREEDOM/pl362/schweizergletschern/")

sf_use_s2(FALSE)

Spherical geometry (s2) switched off

library(ggplot2)

Warning: package 'ggplot2' was built under R version 4.3.2

# Glacier Data Inventories
ch1850 <- st_read(dsn="./data/inventory_sgi1850_r1992/SGI_1850.shp", quiet = TRUE)
ch1931 <- st_read(dsn="./data/inventory_sgi1931_r2022/SGI_1931.shp", quiet = TRUE)
ch1973 <- st_read(dsn="./data/inventory_sgi1973_r1976/SGI_1973.shp", quiet = TRUE)
ch2010 <- st_read(dsn="./data/inventory_sgi2010_r2010/SGI_2010.shp", quiet = TRUE)

#Selecting Specific Parameters within the Shape Files 
ch1850 <- subset(ch1850, select = c(RivLevel0))
ch1931 <- subset(ch1931, select = c(RivLevel0))
ch1973 <- subset(ch1973, select = c(RivLevel0))
ch2010 <- subset(ch2010, select = c(RivLevel0))

# Shape File Switzerland Level 0
switzerland <- st_read(dsn="./data/gadm41_CHE_shp/gadm41_CHE_0.shp", quiet = TRUE)
switzerland<-st_simplify(switzerland, dTolerance = 0.005)

Warning in st_simplify.sfc(st_geometry(x), preserveTopology, dTolerance):
st_simplify does not correctly simplify longitude/latitude data, dTolerance
needs to be in decimal degrees

#Shape File Switzerland Level 3 
switzerland3 <- st_read(dsn="./data/gadm41_CHE_shp/gadm41_CHE_3.shp", quiet = TRUE)
#Simplifying 
switzerland3<-st_simplify(switzerland3, dTolerance = 0.005)

Warning in st_simplify.sfc(st_geometry(x), preserveTopology, dTolerance):
st_simplify does not correctly simplify longitude/latitude data, dTolerance
needs to be in decimal degrees

switzerland3<- subset(switzerland3, select = c(NAME_2, NAME_3))

4.2 Joining Repetitive Data and Simplifying Shape Files

ch1850b<-st_transform(ch1850, st_crs(switzerland3))
ch1850b<-st_simplify(ch1850b, dTolerance = 0.0005)

Warning in st_simplify.sfc(st_geometry(x), preserveTopology, dTolerance):
st_simplify does not correctly simplify longitude/latitude data, dTolerance
needs to be in decimal degrees

ch1850b$glacial_id<-rownames(ch1850b)

ch2010b<-st_transform(ch2010, st_crs(switzerland3))
ch2010b<-st_simplify(ch2010b, dTolerance = 0.0005)

Warning in st_simplify.sfc(st_geometry(x), preserveTopology, dTolerance):
st_simplify does not correctly simplify longitude/latitude data, dTolerance
needs to be in decimal degrees

ch2010b$glacial_id<-rownames(ch2010b)

4.3 Combining Counties with Glaciers

#Combining the Counties with Glaciers and Joining Repetitive Data
library(units)

udunits database from /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/library/units/share/udunits/udunits2.xml

glacial_county_interesec1850<-st_intersection(switzerland3, ch1850b)

although coordinates are longitude/latitude, st_intersection assumes that they
are planar

Warning: attribute variables are assumed to be spatially constant throughout
all geometries

glacial_county_interesec1850$area_glacier_1850<-set_units(st_area(glacial_county_interesec1850), 'km2')

#Combining the Counties with Glaciers and Joining Repetitive Data
library(units)
glacial_county_interesec2010<-st_intersection(switzerland3, ch2010b)

although coordinates are longitude/latitude, st_intersection assumes that they
are planar

Warning: attribute variables are assumed to be spatially constant throughout
all geometries

glacial_county_interesec2010$area<-set_units(st_area(glacial_county_interesec2010), 'km2')

4.4 Glacier Levels in 1850

glacial_county_interesec_1850df<-st_drop_geometry(glacial_county_interesec1850)
glacial_area_by_county1850<-glacial_county_interesec_1850df%>%
  group_by(NAME_3)%>%
  summarize(area_all_glacier_1850_by_county=sum(area_glacier_1850))

Converting the Data to Square Kilometers

#Calculate area county
switzerland3$area_county<-set_units(st_area(switzerland3), 'km2')
swiss_name3<-subset(switzerland3, select = c(NAME_3, area_county))

Merging the datasets by “NAME_3” and isolating the regions that have glaciers inside them

#Left join
merged_1850<-left_join(swiss_name3, glacial_area_by_county1850, by =c("NAME_3" = "NAME_3"))
merged_1850$area_all_glacier_1850_by_county[is.na(merged_1850$area_all_glacier_1850_by_county)]<-0
merged_1850$area_all_glacier_1850_by_county<-as.numeric(merged_1850$area_all_glacier_1850_by_county)

#Plotting 1850 Glacial Regions with Legend 
ggplot()+
  geom_sf(data=merged_1850, aes(fill=area_all_glacier_1850_by_county))+
   labs(title = "Glacial Levels in Swiss Regions for 1850") +  
  scale_fill_viridis_c(option = "turbo")+
  theme_void()

#Plotting the Maps Together 
library(ggpubr)
fig1_1850<-ggplot()+
  geom_sf(data=merged_1850, aes(fill=area_all_glacier_1850_by_county))+
  scale_fill_viridis_c(option = "turbo")+ 
  guides(fill="none")

fig2_1850<-ggplot()+
  geom_sf(data=swiss_name3)+
  geom_sf(data=ch1850b, fill=alpha("Black", 0.9), color=NA)

ggarrange(fig1_1850, fig2_1850, ncol=2)

4.5 Glacier Levels in 2010

glacial_county_interesec_2010df<-st_drop_geometry(glacial_county_interesec2010)  
glacial_area_by_county2010<-glacial_county_interesec_2010df%>%
  group_by(NAME_3)%>%
  summarize(area_glacial=sum(area))

Converting the Data to Square Kilometers

#Calculate area county
switzerland3$area_county<-set_units(st_area(switzerland3), 'km2')
swiss_name3<-subset(switzerland3, select = c(NAME_3, area_county))

Merging the datasets by “NAME_3” and isolating the regions that have glaciers inside them

#Left join
merged_2010<-left_join(swiss_name3, glacial_area_by_county2010, by =c("NAME_3" = "NAME_3"))
merged_2010$area_glacial[is.na(merged_2010$area_glacial)]<-0
merged_2010$area_glacial<-as.numeric(merged_2010$area_glacial)

#Plotting 2010 Glacial Levels in 2010
fig1_2010<-ggplot()+
  geom_sf(data=merged_2010, aes(fill=area_glacial))+
  labs(title = "Glacial Levels in Swiss Regions for 2010") + 
  scale_fill_viridis_c(option = "turbo")+
  theme_void()
fig1_2010

#Plotting the Maps
library(ggpubr)
fig1_2010<-ggplot()+
  geom_sf(data=merged_2010, aes(fill=area_glacial))+
  scale_fill_viridis_c(option = "turbo")+
  guides(fill="none")

fig2_2010<-ggplot()+
  geom_sf(data=swiss_name3)+
  geom_sf(data=ch2010b, fill=alpha("black", 0.9), color=NA)

ggarrange(fig1_2010, fig2_2010, ncol=2)

5 Comparing the Data

Putting the two Maps side by side show a clear presence of Glacial Regression in Switzerland since data gathering began in 1850. The 2010 map shows clear differences especially in the southeast regions of the country with further substantial changes seen in the southwest, the area that boasts Switzerland’s largest mountains as well as the location of the famous Matterhorn.

fig1_1850<-ggplot()+
  geom_sf(data=merged_1850, aes(fill=area_all_glacier_1850_by_county))+
  labs(title = "Glacial Levels 1850") + 
  scale_fill_viridis_c(option = "turbo")+
  theme_void()

fig1_2010<-ggplot()+
  geom_sf(data=merged_2010, aes(fill=area_glacial))+
  labs(title = "Glacial Levels 2010") + 
  scale_fill_viridis_c(option = "turbo")+
  theme_void()

library(ggpubr)
ggarrange(fig1_1850, fig1_2010, common.legend = TRUE, legend="top", ncol = 2)

6 Focusing on a Specific Area

6.1 Mapping Zermatt

Lets Focus on Zermatt, one of Switzerland’s most famous destinations, and home to the famous “Matterhorn” (The Mountain you see on Toblerone Bars)

Zermatt_1850<-subset(glacial_county_interesec1850, NAME_3=="Zermatt")
Zermatt_2010<-subset(glacial_county_interesec2010, NAME_3=="Zermatt")
Aletsch_1850<-subset(glacial_county_interesec1850, NAME_3 == "Fieschertal")
Aletsch_2010<-subset(glacial_county_interesec2010, NAME_3 == "Fieschertal")

geom_sf(data=merged_2010, aes(fill=area_glacial))+

Zermattdiff<-ggplot()+
  geom_sf(data=switzerland, fill = "White")+
  geom_sf(data=Zermatt_1850, fill = "Green")+
  geom_sf(data=Zermatt_2010, fill = "Blue")+
  labs(title = "Glacier Loss in Zermatt") + 
  theme_void()
Zermattdiff

7 Subtracting 2010 Glacier Levels in Comparison to 1850 Levels

By isolating the regions in Switzerland that actually have glaciers in them, as well as merging the glacier levels in 1850 and 2010 allow us to substract the glacial levels in 2010 to the clearly higher ones from centuries ago. The results show the most about of Glacial Regression to be situation in the southwest, home to some of the highest mountains in continental Europe.

together<-merged_2010
together$diff<-merged_2010$area_glacial-merged_1850$area_all_glacier_1850_by_county

together2<-together
together2$diff[together2$diff==0]<-NA
diff_pic<-ggplot()+
  geom_sf(data=together2, aes(fill=diff))+
  labs(title = "Regional Glacial Regression") + 
  scale_fill_viridis_c(option = "turbo", direction = -1)+
  theme_void()
diff_pic